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Abstract 

The results presented here are part of an ongoing 
research program to develop strain rate dependent 
deformation and failure models for the analysis of 
polymer matrix composites subject to high strain rate 
impact loads. State variable constitutive equations 
originally developed for metals have been modified in 
order to model the nonlinear, strain rate dependent 
deformation of polymeric matrix materials. To account 
for the effects of hydrostatic stresses, which are 
significant in polymers, the classical J 2 plasticity theory 
definitions of effective stress and effective plastic strain 
are modified by applying variations of the Drucker- 
Prager yield criterion. To verify the revised 
formulation, the shear and tensile deformation of a 
representative toughened epoxy is analyzed across a 
wide range of strain rates (from quasi- static to high 
strain rates) and the results are compared to 
experimentally obtained values. For the analyzed 
polymers, both the tensile and shear stress- strain curves 
computed using the analytical model correlate well with 
values obtained through experimental tests. The 
polymer constitutive equations are implemented within 
a strength of materials based micromechanics method 
to predict the nonlinear, strain rate dependent 
deformation of polymer matrix composites. In the 
micromechanics, the unit cell is divided up into a 
number of independently analyzed slices, and laminate 
theory is then applied to obtain the effective 
deformation of the unit cell. The composite mechanics 
are verified by analyzing the deformation of a 
representative polymer matrix composite (composed 
using the representative polymer analyzed for the 


correlation of the polymer constitutive equations) for 
several fiber orientation angles across a variety of strain 
rates. The computed values compare favorably to 
experimentally obtained results. 

Introduction 

NASA Glenn Research Center has an ongoing research 
program to develop computational methods for the 
analysis of polymer matrix composites subject to high 
strain rate impact loads. Under these types of loading 
conditions, the deformation of the composite can be 
highly strain rate dependent and nonlinear, which must 
be accounted for within the analytical model. The 
deformation of polymer composites is ordinarily 
assumed to be linear elastic and independent of strain 
rate in transient dynamic finite element codes used for 
impact analysis. 1 

Polymers are known to have a strain rate dependent 
deformation response that is nonlinear above about one 
or two percent strain. Traditionally, viscoelasticity 
models have been used to capture this behavior. 2 
However, there has been an interest in the research 
community in using constitutive equations developed 
for metals, based on plasticity and viscoplasticity 
approaches, to model the nonlinear, strain rate 
dependent behavior of polymers and polymer matrix 
composites. For example, Sun and co-workers 3 ’ 4 
developed a macromechanical, transversely isotropic 
plasticity theory to analyze the nonlinear deformation 
of polymer composites. Bordonaro 5 adapted the 
viscoplasticity theory based on overstress, originally 
developed for metals, to analyze the nonlinear 
deformation of Nylon 66. Pan and co-workers 6,7 and 
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Hsu, et al 8 developed viscoplasticity theories to analyze 
the nonlinear deformation of polymers including the 
effects of hydrostatic stresses. In these cases, the 
effects of hydrostatic stresses on the nonlinear 
deformation, which unlike in metals are significant for 
polymers 9 , were accounted for by applying variations of 
the Drucker-Prager yield criterion 10 to modify the 
definitions of the effective stress and the effective 
plastic strain. A preliminary effort to account for the 
hydrostatic stress effects on the high strain rate 
deformation of polymers was carried out by the authors 
of this report 11 

In this paper, constitutive equations are developed to 
analyze the nonlinear, strain rate dependent 
deformation of polymeric matrix materials in which the 
effects of hydrostatic stresses on the nonlinear 
deformation response are properly accounted for. The 
equations developed in previous work 11 were modified 
in order to more accurately simulate the multiaxial 
stress states found in a composite material. The tensile 
and shear deformation of a representative toughened 
epoxy polymer is characterized and modeled. The 
implementation of the polymer constitutive equations 
within a mechanics of materials based micromechanics 
model incorporating a fiber substructuring approach is 
presented, which permits the analysis of the strain rate 
dependent, nonlinear deformation of polymer matrix 
composites. The tensile deformation of a carbon fiber 
reinforced polymer matrix composite composed using 
the representative polymer examined for the correlation 
of the polymer constitutive equations will be analyzed 
and compared to experimentally obtained results. 

Polymer Constitutive Equations 

To analyze the nonlinear, strain rate dependent 
deformation of the polymer matrix, the Bodner-Partom 
viscoplastic state variable model 12 , which was 
originally developed for metals, is modified. In state 
variable models, there is no defined yield stress; 
inelastic strains are assumed to be present at all stress 
levels. In the “elastic” range of deformation, the 
inelastic strains are assumed to be very small. State 
variables, which evolve with stress and inelastic strain, 
are defined to represent the average effects of the 
deformation mechanisms. 

For this study, temperature effects are neglected, and 
phenomena such as creep, relaxation and high cycle 
fatigue are not accounted for within the equations. The 
nonlinear strain recovery observed in polymers on 
unloading is not currently simulated, and small strain 
theory is assumed to apply. An important point to note 
is that in the original Bodner model as applied to metals 
each of the state variables and material constants had a 


fairly explicit link to specific deformation mechanisms. 
In the application of the applications to the analysis of 
polymers, the relationship of the various constants and 
variables to specific deformation mechanisms is 
somewhat more phenomelogical. Furthermore, all of 
the nonlinearity and strain rate dependence is assumed 
to be due to irreversible deformations. In reality, the 
nonlinear deformation is most likely a combination of 
damage, reversible deformation mechanisms and 
irreversible deformation mechanisms. However, as will 
be shown, for the deformation modes examined in the 
study the model as presented appears to do an adequate 
job of simulating the nonlinear, strain rate dependent 
deformation of the polymers examined. Furthermore, 
the presented model employs a fairly simple 
formulation with a minimum of material constants 
which are fairly easily obtained. 

In the modified Bodner-Partom model used for this 
study, the components of the inelastic strain rate e \. are 

defined as a function of the deviatoric stress 
components Sy, the second invariant of the deviatoric 
stress tensor J 2 and an isotropic state variable Z which 
represents the resistance to molecular flow in the 
polymer (internal stress) in the form 
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where a is a state variable controlling the level of the 
hydrostatic stress effects and D 0 and n are material 
constants. D 0 represents the maximum inelastic strain 
rate, and n controls the rate dependence of the material. 
The elastic components of strain are added to the 
inelastic strain to obtain the total strain. The term a e 
represents the effective stress state in the material, and 
was modified from the original formulation 11 in order to 
account for the effects of hydrostatic stresses in a 
polymeric material. Based on the formulation used by 
Chang and Pan 7 and Hsu, et al 8 , the effective stress is 
defined as follows in this work in order to account for 
the effects of hydrostatic stresses 

o e =J?J^ + a(j kk ( 2 ) 

where is the summation of the normal stress 
components (equal to three times the hydrostatic stress). 
This formulation of the effective stress is based on the 
Drucker-Prager yield criteria 10 . Under pure shear 
loading, the hydrostatic stress is equal to zero and the 
equation reduces to the original formulation 12 , in which 
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the effective stress was set equal to y 3/ 2 for all 
loading conditions. 

The rate of evolution of the internal stress state variable 
Z and the mean stress effect state variable a are defined 
by the equations 

Z = q{Z 1 -Z)e‘ e (3) 

a = q{a l - a)e e (4) 

where q is a material constant representing the 
“hardening” rate, and Z\ and CXi are material constants 
representing the maximum values of Z and a, 
respectively. The initial values of Z and a are defined 
by the material constants Z 0 and Oq. 

The term e\ in Equations 3 and 4 represents the 

effective deviatoric inelastic strain rate, defined as 
follows 



(5) 


where are the components of the inelastic strain rate 

tensor and £ l m is the mean inelastic strain rate. In the 

original Bodner model 12 , the total inelastic strain and 
strain rate are used in the evolution law and are 
assumed to be equal to their deviatoric values. As 
discussed by Li and Pan 6 , since hydrostatic stresses 
contribute to the inelastic strains in polymers, indicating 
volumetric effects are present, the mean inelastic strain 
rate cannot be assumed to be zero, as is the case in the 
inelastic analysis of metals. An important point to note 
is that in the original Bodner model 11 , the inelastic work 
rate was used instead of the effective inelastic strain 
rate in the evolution equation for the internal stress state 
variable. However, for this work the inelastic strain 
rate was deemed easier to work with from both 
computational and characterization points of view, 
particularly in the incorporation of hydrostatic stress 
effects. Since hydrostatic stress effects were not 
considered in the original Bodner model 11 , the 
evolution equation for a is a new addition to the 
formulation. The state variable a is assumed to evolve 
in the same manner as the state variable Z. By using 
this assumption the value of q used in Equation 3 will 
be the same as the value of q used in Equation 4. 


Details on how to obtain the material constants can be 
found in Goldberg, Roberts and Gilat 11 . 

Polymer Equations Correlation Studies 

A series of shear stress-shear strain curves and 
tensile stress-strain curves have been experimentally 
generated for a representative toughened epoxy resin, 
977-2. The shear tests were conducted at strain rates of 
approximately 9xl0' 5 /sec, 2 /sec and 500 /sec. The 
tensile tests were conducted at strain rates of 
approximately 5. 7x1 O' 5 /sec, 1 /sec and 365 /sec. The 
low strain rate tests were conducted using an Instron 
hydraulic testing machine. The high strain rate tests 
were conducted using a split Hopkinson bar. The 
material constants for the constitutive model were 
primarily determined using the shear data. The tensile 
data were used to determine the initial and final values 
of a. The material constants, determined using the 
procedures described in Goldberg, Gilat and Roberts 11 , 
are as follows: E=3500 MPa, v=0.40, D o =lxl0 6 , 
n=0.73, Z 1= 1546 MPa, Z o =380 MPa, q=101, oti=0.201, 
a o =0.223. 

Experimental and predicted shear stress-strain 
curves are shown in Figure 1 and uniaxial tensile stress- 
strain curves are shown in Figure 2 for the 977-2 
material. As can be seen in the figures, for both tensile 
and shear loading the material has a strain rate 
dependent, nonlinear deformation response. For the 
shear test at high strain rates, the sharp increase in 
stress at the beginning of the loading with negligible 
increase in strain is most likely the result of a lack of 
stress equilibrium at the start of loading. The 
oscillations seen in the tensile response at high strain 
rates are most likely due to the specimen geometry 
leading to the stress waves being visible in the 
response. 

The predicted results match the experimental values 
reasonably well for all strain rates for both tensile and 
shear loading. For the medium strain rate, the stresses 
are slightly underpredicted for the initial portion of the 
nonlinear section of the stress-strain curve. This 
discrepancy may be due to the fact that the values of all 
of the material constants were assumed not to vary with 
strain rate. In actuality, this assumption may not be 
completely accurate. While not shown here, as 
discussed in Goldberg, Gilat and Roberts 11 if the 
hydrostatic stress effects are not properly accounted for 
and the material is characterized based on the shear 
results, the predicted stresses under tensile loading 
would be much higher than the experimental values in 
the nonlinear range. Furthermore, if the value of the 
mean stress effect state variable a was kept constant at 
its final value, the stresses in the early part of the 
loading curve would be over predicted. This result is 
significant in that in previous efforts in the literature to 
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account for the effects of mean stresses in the 
deformation response of polymers, the effects of the 
hydrostatic stresses have been assumed to be constant 
over the entire range of loading. 

Composite Micromechanical Modeling 

Micromechanical techniques are used to predict the 
effective properties and deformation response of the 
individual plies in a composite laminate. The effective 
properties and deformation response are computed 
based on the properties of the individual constituents. 
Lamination theory can then be used to compute the 
effective deformation response of the entire composite. 
The constitutive equations described above have been 
implemented within a micromechanics method in order 
to enable the prediction of the nonlinear, strain rate 
dependent deformation response of polymer matrix 
composites with the effects of hydrostatic stresses 
incorporated into the analysis. The micromechanics 
method has been described extensively in Goldberg 13 . 

A summary of the methodology will be given here. 

For this work the unit cell, the smallest material unit 
for which the response can be considered to be 
representative of the entire composite ply, is defined to 
consist of a single fiber and its surrounding matrix. 

Due to symmetry, only one-quarter of the unit cell was 
analyzed. The composites are assumed to have a 
periodic, square fiber packing and a perfect interfacial 
bond is specified. The fibers are assumed to be 
transversely isotropic and linear elastic with a circular 
cross-section. The matrix is assumed to be isotropic, 
with a rate dependent, nonlinear deformation response 
computed using the equations described in the previous 
section of this report. A key assumption of this 
approach is that the in-situ matrix properties are 
equivalent to the bulk properties of the polymer. 
However, the advantage of using this type of 
methodology is that it is simpler to conduct experiments 
on the pure resin and to determine the material 
constants from the pure resin data as opposed to trying 
to back out the resin properties from composite test 
data. Furthermore, a key goal of this research is to 
provide a methodology that facilitates reducing the 
amount of testing of the composite that is required to 
obtain strain rate dependent material properties that can 
be input into a finite element code. Conducting strain 
rate dependent tensile tests on the pure resin and using 
that data to predict the composite deformation response 
is also much simpler than conducting tests on the 
composite. However, if in comparing test data obtained 
from composite specimens to analytical predictions it 
appears that the bulk matrix properties do not 
accurately reflect the in-situ state of the matrix, the 
polymer properties can always be appropriately 
adjusted. 


The unit cell is divided up into an arbitrary number 
of rectangular, horizontal slices of equal thickness, as is 
shown in Figure 3. Similar approaches have been used 
by researchers such as Whitney 14 , Greszczuk 15 and 
Mital, et al 16 . Each slice is assumed to be in a state of 
plane stress. This assumption is made based on the fact 
that laminate theory will be applied to each ply of the 
composite laminate, which implies that the unit cell and 
every slice within the unit cell is in a state of plane 
stress. The top and bottom slices in the unit cell are 
composed of pure matrix. The remaining slices are 
composed of two subslices; one composed of fiber 
material and the other composed of matrix material. 

For the slices containing both fiber and matrix, the out- 
of-plane stresses can be nonzero in individual subslices, 
but the volume average of the out-of-plane stresses 
must be equal to zero. By using this approach, the 
behavior of each slice is decoupled, and the response of 
each slice can be determined independently, which 
significantly reduces the level of complexity in the 
analysis. Laminate theory is then used to obtain the 
effective response of the unit cell. 

The thickness, fiber volume ratio and thickness ratio 
(the ratio of the slice thickness to the total unit cell 
thickness) for each slice can be determined using the 
composite fiber volume ratio and geometric principles. 
The unit cell is assumed to measure one unit in length 
by one unit in height. The first step is to compute the 
area of the cross-section of the fiber within each slice. 
The overall diameter of the fiber (df) is related to the 
fiber volume fraction of the overall composite (V f ) 
through the following relationship 



and this term can be used along with the standard 
geometric definition of the radius of a circle to compute 
the horizontal coordinate of any point on the outer 
surface of the fiber in terms of the fiber volume fraction 
and the vertical coordinate. The area of the portion of 
the fiber contained within each slice (A/) within the 
one-quarter of the unit cell which is analyzed can 
computed by integrating the resulting expression 
between the vertical (z) coordinates of the top and 
bottom of slice “i” 



which is also the equivalent area of the rectangular fiber 
slice in the portion of the unit cell which is analyzed 
(one-quarter due to symmetry). 
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The fiber volume fraction of each slice composed of 
fiber and matrix is equal to the fiber area in each slice 
divided by the total slice area. The thickness ratio for 
each slice composed of both fiber and matrix is equal to 
the slice thickness divided by the assumed total height 
of the analysis cell. The fiber volume fraction of the 
top slice consisting of matrix only is equal to zero, and 
the thickness ratio of the top slice is equal to one minus 
the sum of the thickness ratio of the remaining slices. 

The effective properties, effective inelastic strains 
and effective thermal strains of each slice are computed 
independently. The responses of each slice are 
combined using laminate theory to obtain the effective 
response of the corresponding lamina. Micromechanics 
equations are developed for those slices composed of 
both fiber and matrix material. The stresses in the 
slices composed of pure matrix can be computed using 
the matrix elastic properties and inelastic constitutive 
equations. The standard transversely isotropic 
compliance matrix (or isotropic in the case of the 
matrix) is used to relate the local strains to the local 
stresses in the fiber and matrix. Each slice is assumed 
to be in a state of plane stress on the global level, but 
out-of-plane normal stresses can exist in each subslice. 
Along the fiber direction, the strains are assumed to be 
uniform in each subslice, and the stresses are combined 
using volume averaging. The in-plane transverse 
normal stresses and the in-plane shear stresses are 
assumed to be uniform in each subslice, and the strains 
are combined using volume averaging. The out-of- 
plane strains are assumed to be uniform in each 
subslice. The volume average of the out-of-plane 
stresses in each subslice is assumed to be equal to zero, 
enforcing a plane stress condition on the global level 
for the slice. 

An orthotropic compliance matrix is used to relate 
the strains (By) to the stresses (cTy) in each constituent. 
The addition of the inelastic strain components to the 
standard orthotropic elastic constitutive law facilitates 
the incorporation of inelasticity into the constitutive 
relations. For the fiber, which is assumed to be linear 
elastic, these components are neglected. 

By combining the uniform stress and uniform strain 
assumptions with the constituent stress-strain relations, 
a system of four simultaneous equations results that can 
be solved for the unknown stresses in the subslices. 

The total strains and subslice inelastic strains are 
considered to be the known values in solving this 
problem. By substituting the subslice stresses back into 
the equations defining the uniform stress assumptions, 
the effective elastic constants, effective inelastic strains 
and effective thermal strains of the slice can be 
computed. By applying classical laminate theory at this 
point, the effective stiffness matrix, effective inelastic 
strains and effective thermal strains for the unit cell are 
computed. Laminate theory is applied once again to 


obtain the effective properties and force resultants due 
to inelastic and thermal strains for the multilayered 
composite laminate. Further information on all of these 
procedures can be found in Goldberg 13 . 

Simulation of Strain Rate Dependent Composite 
Deformation 

To verify the micromechanics equations and the 
implementation of the polymer constitutive equations 
within them, a series of analyses have been carried out 
using a representative polymer matrix composite 
system that exhibits a strain rate dependent, nonlinear 
deformation response. The material examined consists 
of carbon IM7 fibers in the 977-2 toughened epoxy 
matrix discussed earlier. Longitudinal tensile tests were 
conducted on composite laminates with various fiber 
orientations. Tests were conducted at strain rates of 
about 5xl0' 5 /sec, about 1.0 /sec and about 400-600 
/sec. Dog-bone shaped specimens were used with a 
gage length of approximately 0.9525 cm. The low 
strain rate testing was conducted using an Instron 
hydraulic testing machine. The high strain rate tests 
were conducted using a tensile split Hopkinson bar 
apparatus. 

The IM7/977-2 composite has a fiber volume ratio 
of 0.60. The material properties used in this study for 
the IM7 fiber include a longitudinal modulus of 276 
GPa, a transverse modulus of 13.8 GPa, a longitudinal 
Poisson’s ratio of 0.25, a transverse Poisson’s ratio of 
0.25 and an in-plane shear modulus of 20.0 GPa. These 
properties are as given in Gates, et al 17 , with the 
exception of the value for the transverse Poisson’s ratio 
was taken from Murthy, et al 18 based on representative 
carbon fiber data. The material properties of the 977-2 
resin are as given before. 

Experimental and computed longitudinal tensile 
stress-strain curves for two laminate configurations 
([45°] and [±45°] s ) of the IM7/977-2 material are 
shown in Figure 4 and Figure 5. These laminate 
configurations were chosen due the pronounced 
nonlinearity and strain rate dependence observed in the 
experimental results. Five fiber slices were used in the 
portion of the unit cell which was analyzed. This value 
was found to yield sufficiently converged answers. In 
Figure 4, results for the [45°] laminates at strain rates of 
4.75xl0' 5 /sec, 1.2 /sec and 405 /sec are shown. In 
Figure 5, results for the [±45°] s laminates at strain rates 
of 9xl0* 5 /sec, 2.1 /sec and 604 /sec are shown. 

As can be seen in the figures, the analytical model 
captures the strain rate dependence and nonlinearity 
observed in the experimental stress-strain curves. For 
the low and moderate strain rate curves for both 
laminates and the high strain rate curve for the [45°] 
laminate, the comparison between the experimental and 
analytical results is reasonably good. The stresses at 
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the moderate strain rate for the [±45°] s laminate are 
somewhat overpredicted and the stresses at the low 
strain rate for the [45°] laminate are somewhat 
underpredicted for reasons which are clear at this time. 
For the [±45°] s laminate at high strain rates, the stresses 
are significantly overpredicted, indicating that there 
may have been some damage in the experimental 
specimens which is not predicted with the current 
analytical model. Alternatively, there could be an error 
in the material constants. The specific causes for the 
discrepancies need to be investigated further. 

Conclusions 

A set of constitutive equations have been developed 
to analyze the strain rate dependent, nonlinear 
deformation of polymeric matrix materials, including 
hydrostatic stress effects. The tensile and shear 
deformation response of a representative polymer over 
a range of strain rates have been successfully predicted. 
The constitutive equations have been implemented 
within a strength of materials based micromechanics 
approach in which the unit cell is subdivided into a set 
of independently analyzed slices. The micromechanics 
technique was used to predict the strain rate dependent 
deformation of a representative polymer matrix 
composite for two fiber orientations. Overall, the 
experimental stress-strain curves were accurately 
predicted using the model, with some discrepancies for 
the high strain rate results which need to be investigated 
further. 
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Figure 1: Experimental and computed shear stress-shear strain curves for 977-2 resin at strain rates of 9xl0' 5 /sec 
(Low Rate), 1.9 /sec (Medium Rate) and 518/sec (High Rate). 
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Figure 2: Experimental and computed tensile stress-strain curves for 977-2 resin at strain rates of 5.7xl0' 5 /sec (Low 
Rate), 1.9 /sec (Medium Rate) and 365 /sec (High Rate). 
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Figure 3: Schematic showing relationship between unit cell and slices for micromechanics. 
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Figure 4: Experimental and predicted results for IM7/977-2 [45°] laminates at strain rates of 4.7xl0' 5 /sec (Low 
Rate), 1.2 /sec (Medium Rate) and 405 /sec (High Rate). 
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Figure 5: Experimental and predicted results for IM7/977-2 [±45°] s laminates at strain rates of 9xl0' 5 /sec (Low 
Rate), 2.1 /sec (Moderate Rate) and 604 /sec (High Rate). 
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